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Abstract 

We give an review the HBT puzzle, and argue that its resolution requires the introduction of new 
physics close to the phase transition scale. We argue that a candidate for this new physics is bulk 
viscosity, recently postulated to peak, and even diverge, close to the phase transition temperature. 
We show that such a viscosity peak can force the system created in heavy ion collisions to become 
unstable, and filament into fragments whose size is weakly dependent on the global size of the 
\Jq ■ system, thereby triggering freeze-out. 

<N " 

1 The HBT puzzle 

o, 

ON One of the most unexpected, and as yet unexplained, experimental results found at the Relativistic 
^ Heavy Ion Collider (RHIC) concerns the description of particle interferometry observables [lj. It was 
\ originally expected that the deconfined matter would be a highly viscous, weakly interacting quark 
^ gluon plasma [2]. Thus, ideal hydrodynamics would not provide a good description of flow observables 
o3 sensitive to the early stages of the collision, such as azimuthal anisotropy. The signature of choice of 
a phase transition from hydrodynamics models, one less sensitive to viscosity, was to be an increase of 
the "out" to "side" emission radius ratio (R Q and see Fig. [1] center panel) due to longer lifetime of 
the system, caused by the softening of the equation of state in the transition/crossover region [3j. 

The data, however, exhibited an opposite behaviour. Hydrodynamic simulations provided a good 
description of transverse momentum spectra and their azimuthal anisotropy. The same simulations, 
however, failed to describe HBT data Perhaps the most surpising aspect of the problem is the 
way in which the data does not fit: Measured parameters R Q and R s are nearly identical. (Fig. [I] left 
and right panel) Their (positive) difference R 2 Q — R 2 S is thought to depend on the duration of particle 
emission. Hence, it looks like the fireball emits particles almost instantaneously and does not show any 
sign of phase transition or crossover. Hydrodynamics, with "reasonable" freeze-out condition (such as 
a freeze-out temperature of 100 MeV or so) can not describe this even qualitatively. This behaviour, 
when compared to lower energy data, exhibits remarkably good scaling with multiplicity (Fig. [1] left 
panel). The scaling's existance, however, is by itself surprising since the QCD equation of state, with 
it's critical density for a phase transition, should break it. 
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Figure 1: Dependence of HBT radii on multiplicity, across different energies [lj (left panel), an ex- 
pectation from hydrodynamics [3] (center panel), and a fit with hydrodynamic model [4j (right panel) 



There could be three possible approaches to the HBT puzzle. It could be that the system is simply 
too complicated, and that once we include all possible improvements (full 3D calculation, viscosity, 
hadronic kinetic afterburner, in-medium hadron modifications, pre-existing flow etc.), everything will 
fit. It could also be that we are drastically misunderstanding the data, and the HBT puzzle is a 
symptom of inapplicability of hydrodynamics to heavy ion collisions. Finally, it could be that the 
hydrodynamic approach is basically correct, but just one element of relevant physics is omitted [5] 

The second possibility is unlikely because, in some ways, hydrodynamic prescription does fit HBT 
data. Scaling of HBT radii with the multiplicity rapidity density (dN ' j 'dy) 1 ^ , over a large range of 
energies [1] is typical for an isentropically expanding fluid that suddenly breaks apart. In addition, the 
good description, within parameters compatible with what is needed to describe flow, of the azimuthal 
dependence of HBT radii [4J, also suggests that the hydrodynamic framework is a good ansatz for 
describing the matter produced in heavy ion collisions up to freeze- out. 

The first possibility also appears problematic: successful models and/or parametrisations of the 
freeze-out which describe HBT radii are found in the literature [6], and they could provide a way to 
gain insights into what is missing. However, we feel that successful description involves a dynamical 
description from initial conditions plus a freeze- out criterion, rather than a fit to data with assumptions 
put in "by hand" . Such a description is so far lacking. Furthermore, the most plausible refinements to 
hydrodynamics, namely implementation of fully three-dimensional models and the addition of a kinetic 
theory afterburner [7] do not do anything to solve the HBT discrepancy, suggesting that the problem 
is not refinements but rather one large missed physical effect 



2 Bulk viscosity close to T c 

The bulk viscosity of high temperature strongly interacting matter has recently been calculated using 
perturbative QCD [8], and found to be negligible, both in comparison to shear viscosity and w.r.t. 
its effect on any reasonable collective evolution of the system. This is not surprising: The QCD 
Lagrangian, as long as no "heavy" quarks are present, is nearly conformally invariant [8]. Since, within 
a fluid, the violation of conformal symmetry is linearly proportional to a bulk viscosity term, the near 
conformal invariance of the QCD Lagrangian should guarantee that bulk viscosity is nearly zero, in 
the perturbative regime. In the hadron gas phase, the numerous scales associated with hadrons render 
conformal invariance a bad symmetry, and hence it is natural that bulk viscosity not be negligible. 
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This is, again, rooted in a fundamental feature of QCD: the non-perturbative conformal anomaly, 
that manifests itself in the scale (usually called Kqcd) at which the QCD coupling constant stops being 
small enough for the perturbative expansion to make sense. This scale coincides with the scale at which 
confining forces hold hadrons together. This violation of conformal invariance is not seen perturbatively, 
but should dominate over the perturbatively calculated bulk viscosity as temperature drops close enough 
to the QCD phase transition. What happens to bulk viscosity in this regime, where hadrons are not yet 
formed, presumably the matter is still deconfined, but conformal symmetry is badly broken? There is 
a quenched lattice calculation, and compelling physical arguments [9] that bulk viscosity rises sharply, 
or even diverges, close to the phase transition temperature. 

Lattice simulations find that Tff (=0 for a conformally invariant system), increases rapidly close to 
T c . Remembering that the shear (77) and bulk (£) viscosities roughly scale as [10] 

V ~ ^elastic? 14 C ~ Q - vfj Ti ne i as ti c r 4 (l) 

where T(i ne ) e i ast ic refers to the equilibration timescale of (ine) elastic collisions. Finite temperature sum 
rules in conjunction with lattice data [9j give a sharp rise in bulk viscosity. 

The rise is, in fact, likely to be considerably sharper than [9] suggests. The dependence of r ine i ast i C on 
temperature can be guessed from the fact that, at T c , the quark condensate (qq) acquires a finite value, 
and the gluon condensate (G^G^) sharply increases at the phase transition. "Kinetically" , therefore, 
timescales of processes that create extra qq and GG pairs should diverge close to the phase transition 
temperature, by analogy with the divergence of the spin correlation length in the Ising model close to 
the phase transition. The sharp rise of bulk viscosity can also be understood within string kinetics: 
confinement, microscopically, can be thought of as a "string tension" appearing in the potential. Even 
in a regime where the momentum exchange of the average collision is more than enough to break the 
string, and the relevant degrees of freedom are still quarks, not mesons, string tension introduces a huge 
change in kinetics: what were elastic collision without string tension, become inelastic onece string 
tension appears. Even if the energy needed to break the string is low, over many collisions, the heat 
energy would be converted into creating more slightly colder, less pressing particles. That's exactly the 
kind of processes that contribute most to bulk viscosity [TT] . 

These arguments give evidence to the conjecture that, close (from above) to T c , bulk viscosity goes 
rapidly from a negligible value to a value capable of dominating the collective evolution of the system. 
For our analysis, we assume the bulk viscosity to be of the form ( = s (z p qcd + ex P( — ^2 )) where 
t = T — T c and a = 0.01T C and z p qcd ~ 10 -3 [8j. At T > T c , this ansatz provides a reasonable fit to 
the results of [9], considering the peak height and width of the distribution are still unknown. 



3 Bulk viscosity-triggered fragmenting 

To study the effect of our conjectured behaviour of bulk viscosity on solutions to hydrodynamics, we 
perform a stability analysis [12] of a boost-invariant solution to the Navier-Stokes equations (both the 
ID and 3D cases). The Navier-Stokes equations with boost-invariant symmetry [21 [[21 US] can be 
rewritten in terms of the Reynolds number i?, the entropy s the co-moving time r, the total number 
(N) of dimensions, and the dimensionality of the homogeneous expansion (M) 

,- M d{r M s) Ms 

M = 1 N = 3 corresponds to Bjorken hydrodynamics [T3j , M = N = 3 to a homogeneus 3D expansion. 
The Reynolds number is a function of temperature T, bulk and shear viscosity C and rj and entropy 
density s: R' 1 = (2(1 - M/N)r) + MC)/(Tst). With given expressions for s(T), r)(T), ((T), this set of 
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equations becomes closed and solvable. For the equation of state, we use a parameterization of lattice 
data (We have checked that our results do not vary qualitatively if the ideal EoS is used) . 

We follow the stability analysis performed in [12]. The amplitude of a generic perturbation to the 
ID Boost-invariant system is a vector x in the two dimensional space of entropy perturbations and flow 
(rapidity y) perturbations, and its frequency in rapidity can be decomposed into Fourier components 

5s 

£(y) = %( k ) elky : Xi = — , X 2 = y~ Vspacetime (3) 
k S 



The equation of motion for x will then be given by 

d ( x, \ = / A n A 12 
\ A 2 i A 22 





(4) 



where A^ is a real matrix function of s, R and k. We refer the reader to Eq. 4.23-4.26 of p2] for the full 
form of Aij. The system's stability can be understood via the behavior of the modulus of x, X = oFx 

Q — » — » — > —>T 

t—X = x T M tJ x , Mjj = Aij + A %3 (5) 

Since is real and symmetric, it will always have two real eigenvalues, X max and X min (corresponding 

eigenvectors x max ^ min ), as well as orthogonal matrices B^ diagonalizing it. Defining $ = Xjwe see 
that X min y min < < ^rnaxVmax Thus, if X min > 0, the system is unstable, since perturbation in any 
direction will produce a positive growth rate. An instability will grow as a power-law with t/t^, where 
Ti n i is the starting time of the perturbation x{r) ~ (r ^/t^) a x(t^). If \ ma x < 0, on the other hand, the 
system will be stable against perturbations, an instability will be suppressed as ~ (t / T ini ) x ( <0 \ 

If Xmin < and X max > some modes will be stable and some will be unstable. In the latter case, the 

time dependence of A^ will in general continuously rotate x min and x max in time, stopping the growth 

of the instability: Since A^ is time-dependant, x max (T in i) might be in the direction of x min (r > r ini ) a 
short time later. Solving Eq. [Hwill take this effect into account. 

In what follows, we use z = 0, 0.1, 1 T c and evolve the system from an initial temperature T = 0.3 
GeV and comoving time 0.6 fm. Note that the qualitative features of our study are independent of the 
details of the evolution before T c , such as the initial temperature and timescale. 

Figure [2] (upper panel) shows the temperature, entropy density and total entropy in the central 
rapidity unit as a function of time. As soon as z becomes non-negligible (i.e. viscous forces dominate 
around T c ), the kinematic evolution of the system "freezes". The system then stays at nearly constant 
temperature, through it keeps producing more and more entropy at the expense of advective energy. 
At first sight, large values of z are excluded by HBT data and multiplicity measurements. However, 
we will show that this long phase is unstable against small perturbations. Thus, its further evolution 
will not be given by the background solution, but by a rapid formation of local inhomogeneities. 

Figure [2] (second panel) shows the X min and X max Eigenvalues corresponding to representative k = 2, 8 
(other values of k were checked not to vary significantly wrt those presented here). As can be seen, 
the peak in bulk viscosity forces the growing/damping rates to increase rapidly. Thus, any initial 
perturbation in the unstable direction will rapidly grow to a value comparable with the background, 
unless the system's evolution will stop the growth by rotating the direction of the unstable modes. 

Figure [2] (third panel) examines whether this occurs for larger values of z . If the peak of viscosity 
is negligible, the unstable eigenvector keeps rotating throughout the evolution of the system. Thus, 
even an unstable mode's growth will very quickly stop growing since the dynamics will turn it into 
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Figure 2: Instability analysis of the boost-invariant solution. Top panel shows the evolution of the 
background temperature (solid), entropy density (dot-dashed) and total entropy (dashed). The two 
middle panels show the eigenvalues and eigenvectors of the k = 2 (solid) and k = 8 unstable mode. 
The bottom panel shows the full evolution of the linearized instability equation of motion. The three 
columns represent various appearances of the viscosity peak. See text for details. 



a damped mode. When z dominates, however, when T approaches T c , the direction of the unstable 
modes experiences an abrupt rotation. Then it stays constant throughout the time the system travels 
through the viscosity peak (this time increases strongly as z increases), and gets rotated again as the 
peak is passed. The reason for this behavior is clear from the "freeze" of the background evolution in 
the top panel. Thus, at large z unstable modes have all the required time to grow. 

Finally, Fig. [2] (bottom panel) shows the explicit solution of Eq. ES At each time-step, a perturbation 
is born in the unstable eigenvector mode, and then evolved until the end of the evolution of the system. 
The plot shows Z = X (t) / X (T ini ) , the ratio of strength of the perturbation to the initial strength as a 
function of time (Note that the starting value is always unity). A large Z does not mean the evolution 
equation is invalid: For any point in the graph, there will be a small enough perturbation amplitude 
that survives as a perturbation in the subsequent evolution. The probability of a larger perturbation 
forming and significantly modifying the background, however, should grow strongly with Zq. It is clear 
that any microscopic mechanism seeding instabilities at the scale Z ~ 10 _1 , uniformly distributed in x 
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and at a rate of ~ fm~ 4 is likely to generate power-law growing instabilities a few fm after T ~ T c . 
These instabilities should reach Z > 1, and hence play an important role in the subsequent evolution 
of the system, a few fm after that. 

The resulting scenario is phenomenologically similar to the fragmenting induced by supercooling in 
a first order transition [14], but originating from hydrodynamic rather than thermodyamic instabilities. 

We do not expect second order hydrodynamics p3] to play a big role in starting the instabilities: 
As argued in [5j, as long as the system's expansion rate as the system approaches T c is smaller than the 
relaxation time m, r n^^ < 1 second order terms will not prevent viscous corrections to the pressure 
of the order of £/r, but merely localize their propagation. The (admittedly unreliable) estimates from 
strongly coupled CFT [16] suggest that this criterion is amply satisfied, especially considering that the 
increase in viscosity causes the background solution to slow down over a timescale much bigger than 
Tn- Second order hydrodynamics might however play a dominant role in stopping the instabilities once 
they grow to a value comparable to the background. 

Recent numerical simulations, with Israel-Stewart hydrodynamics, provide evidence that this de- 
scription is true, as cluster-like instabilities seem to be present (Fig. 11 of [17] ) 

4 Fragmenting and the HBT puzzle 

HBT radii are directly related to the system's spacetime correlation tensor [H [lf*l 



while R 2 S (K) is simply given by ((Ax s ) 2 ). Here the K is the sum of the two four-momenta of the pion 
pair. As remarked in [I], the R Q ~ R s result is not easy to reconcile with naive hydrodynamics plus a 
straight-forward (critical temperature) emission for several reasons; 

Firstly, the higher the energy, the longer the emission time, the larger is the expected discrepancy 
between R Q and R s . If the system starts close to the mixed phase, the timescale of freezing out should 
be longer still due to the softest point in the equation of state. Hence, a generic prediction from Eq. 
([6]) is that R /R s > 1, largely increases with energy, and exhibits a peak when the energy density is 
such that the system starts within, or slightly above the mixed phase. This is in direct contrast with 
experimental data, where R /R s — 1 is a feature at all reaction energies. 

Secondly, generally (Ax At) < 0, since particles on the outer layer freeze-out first. This increases 
R /R s further (cf. eq[6]). Time dilation due to transverse flow does reverse this dependence [I]. 

Fragmentation of the bulk could help solving these problems. Firstly, fragment size, density and 
decay timescale, is approximately independent of either reaction energy or centrality. Hence, the near 
energy independence of the (comparatively short) emission timescale, and hence of R Q /R S) should be 
recovered. Secondly, if the decay products do not interact (or do not interact much) after fragment 
decay, it can also be seen that (Ax At) can indeed be positive: outward fragments are moving faster, 
resulting in time dilation. This effect can be offset by time dilation of fragment decay by increasing the 
temperature at which fragments form, or by increasing fragment size. Recovering the linear scaling of 
the radii with (dN/dy) 1 / 3 ^ 7Vf ragments ) [I], while maintaining the correct R /R s is also possible if the 
fragments decay when their distance w.r.t. each other is still comparable to their intrinsic size. 

The bulk- viscosity-driven freeze-out adds another parameter to ab initio HBT calculations: in ad- 
dition to critical temperature/energy density, we now have the fragment size. To see whether this 

1 Here I ("long") is the z direction (parallel to the beam), o ("out") is the direction of the pair momentum, and s 
("side") is the cross product of the previous two. Averaging is done over the emission function S(x,p), the probability of 
producing a particle of momentum p at x. See [1 for details 
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helps solving the HBT problem, output from hydrodynamics with a high (T ~ T c ) freeze-out tempera- 
ture should be fragmented into piecess with a certain distribution in size, which then produce hadrons 
according to the prescription in the Appendix of [5]. 

In conclusion, we have introduced the HBT puzzle, and explained how bulk viscosity could help 
solve it by triggering the fragmenting of the system close to the critical temperature. We hope that, in 
the near future, this scenario will be developed to the point of experimental falsification. 
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